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Abstract 

We present the two-loop virtual QCD corrections to the production of heavy quarks in gluon fusion. 
The results are exact in the limit when all kinematical invariants are large compared to the mass of 
the heavy quark up to terms suppressed by powers of the heavy-quark mass. Our derivation uses 
a simple relation between massless and massive QCD scattering amplitudes as well as a direct 
calculation of the massive amplitude at two loops. The results presented here together with those 
obtained previously for quark-quark scattering form important parts of the next-to-next-to-leading 
order QCD corrections to heavy-quark production in hadron-hadron collisions. 



1 Introduction 



The production of heavy quarks at hadron colliders is an important process for a range of reasons. 
Experimentally, events with a top-quark, being the heaviest quark known thus far, lead to very 
characteristic signatures in a hadronic scattering process. This allows for event reconstruction in 
a variety of channels, which e.g. at LHC is supplemented by an anticipated very high statistics 
for the production cross section [1,2]. In this way precision measurements of the top-mass, of 
(differential) ?f-distributions and also tests of the production and the subsequent decay mechanism 
(including anomalous couplings and top-spin correlations) can be conducted at Tevatron or LHC. 
The precise knowledge of the top-quark parameters has direct impact on precision tests of the 
Standard Model, in particular on the Higgs sector due to the large top- Yukawa coupling. Moreover, 
events with top-quarks will also make up for a large part of the background in searches for the 
Higgs boson or new physics. 

In the case of bottom-quark production important measurements are inclusive 5-meson and b- 
jet production at moderate to large transverse momentum. At colliders Z^-flavored jets are produced 
for instance in decays of top-quarks, the Higgs boson and numerous particles proposed in exten- 
sions of the Standard Model. Differential b-jet distributions are also important for measurements of 
parton distributions. Predictions for heavy-quark hadro -production have theoretical uncertainties, 
which can even be larger than e.g. the corresponding predictions for the production of light jets. 
This is due to radiative corrections in Quantum Chromodynamics (QCD), where typically large 
logarithms at higher orders appear. At present, the next-to-leading order (NLO) QCD corrections 
for heavy-quark hadro-production are known [3-9] and, once they become available, the com- 
plete next-to-next-to-leading order (NNLO) QCD corrections will reduce the uncertainty of theory 
predictions by improving the stability with respect to scale variations. For inclusive 5-meson pro- 
duction [10, 11], for instance, the formalism of perturbative fragmentation functions of a heavy 
quark [12] has already been extended through NNLO by providing the initial conditions [13, 14] 
together with the time-like (non-singlet) evolution [15]. Alternative ways devoted to control the 
theoretical uncertainty in the case of Z?-quarks, for example through the definition of dedicated jet 
algorithms have been proposed recently [16]. 

In the hadronic production of heavy quarks large perturbative QCD corrections arise in dif- 
ferent kinematical regions. If the heavy quarks are produced with (partonic) center-of-mass en- 
ergies s close to threshold, then Sudakov-type logarithms, powers of log(|3), appear to all orders 
in perturbation theory, where |3 = \/\ — Am~ / s is the velocity of the heavy quark. The necessary 
resummation has been carried out to next-to-leading logarithmic (NLL) accuracy [17-19] and has 
successfully improved the phenomenology of top-quark production at Tevatron. In the high-energy 
regime on the other hand, when the heavy quarks are fast, logarithms in the quark mass m appear. 
These log(7M)-terms are of coUinear origin and dominate cross sections when m becomes negli- 
gible in comparison to other kinematical invariants, such as e.g. for the inclusive Z?-jet spectrum 
at large transverse momentum. All-order predictions can be achieved by means of an explicit re- 
summation. The necessary technology for resumming both the incoming and outgoing coUinear 
logarithms to NLL accuracy is available, see e.g. Refs. [20,21] and references therein. 
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In this article, we want to further improve the precision of fixed-order perturbative QCD cor- 
rections. To that end, we present results for the virtual QCD corrections at two loops for the 
pair-production of heavy quarks in gluon fusion. Together with the corresponding results for 
quark-quark scattering obtained previously by us [22] these are essential parts of the complete 
NNLO QCD corrections. To be precise, we calculate the interference of the two-loop with the 
Bom amplitude. We work in the limit of fixed scattering angle and high energy, where all kine- 
matic invariants are large compared to the heavy-quark mass m. Thus, our result contains all 
logarithms log(m) as well as all constant contributions (i.e. the mass-independent terms) and we 
consistently neglect power corrections in m. 

In our calculation we employ two different methods. On the one hand, we apply a general- 
ization of the infrared factorization formula for massless QCD amplitudes [23,24] to the case of 
massive partons [25]. In a nut-shell this results in an extremely simple universal multiplicative 
relation between a massive QCD amplitude in the small-mass limit and its massless version [25]. 
In this way, we can largely use for our derivation the results of the NNLO QCD corrections to 
massless quark-gluon scattering (i.e. gg qq). These have been computed for the squared ma- 
trix element, i.e. the interference of the two-loop virtual corrections with the Bom amplitude in 
Ref. [26], while results for the individual (independent) helicity configurations of the two-loop 
amplitude itself have been given in Refs. [27, 28]. On the other hand, we perform a di- 

rect calculation of the relevant Feynman diagrams in the massive case followed by a subsequent 
expansion in the small-mass limit by means of Mellin-Barnes representations. 

The outline of the article is as follows. In Section[2lwe give some basic definitions and present 
a short summary of our methods in Section [3l There we briefly explain the essence of the QCD 
factorization approach to calculate massive amplitudes and give the relevant formulae. We also 
highlight the key steps of the direct evaluation of Feynman diagrams in the small-mass limit and, in 
particular, comment on the non-planar topologies. Section |4] contains our results and we conclude 
in Section [51 Appendix [Al gives the explicit result for a massive non-planar scalar integral in the 
small mass limit. 



2 Setting the stage 

The pair-production of heavy quarks in the gluon fusion process corresponds to the scattering 

8{pi)+g{P2) Q{P3,in) + Q{p4,m), (1) 

where pi denote the gluon and quark momenta and m the mass of the heavy quark. Energy- 
momentum conservation implies 

Pt + P'2=P3+P4- (2) 

Following the notation of Ref. [26] we consider the scattering amplitude for the process ([T]) at 
fixed values of the external parton momenta pi, thus pl= P2 = and p^ — p^ — m^- The amplitude 
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M may be written as a series expansion in the strong coupling tts, 

|f^(o)) + (|)|^0)) + (|)V(^^) + o(a3) 



Ana. 



(3) 



and we define the expansion coefficients in powers of CXs(aj^) / (27r) with jj. being the renormalization 
scale. We work in conventional dimensional regularization, J = 4 — 2£, in the MS-scheme for the 
coupling constant renormalization. The heavy mass m on the other hand is always taken to be the 
pole mass. 

We explicitly relate the bare (unrenormalized) coupling to the renormalized coupling tts by 



po/as\ , (^l Ipi 



In 



2 £ 



In 



(4) 



where we put the factor 5e = (47l)^exp(— £'^) = 1 for simplicity and (3 is the QCD p-function 
[29,30] 



(5) 
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The color factors are in a non-Abelian SU(//) -gauge theory =N,Cp = {N^ —l)/2N and Tp = 
1/2. Throughout this letter, A'^ denotes the number of colors and Hj- the total number of flavors, 
which is the sum of Hj light and heavy quarks. 

In the following, we will confine ourselves to the discussion of the squared amplitude for the 
process ©, although it should be clear that our approach and the results of the present article can be 
easily extended to the (color ordered) partial amplitudes for the individual helicity combinations of 
the massive two-loop amplitude |5Vf (2^) itself. There one would rely in particular on Refs. [27,28]. 

For convenience, we define the function (£, m, sJ./j) for the squared amplitudes summed over 
spins and colors as 



£|fAf(g + g^<2 + <2)P = ^(£,m,5,?,A/). 
^ is a function of the Mandelstam variables s, t and u given by 

s = {p\+P2f, t = {pi-p2,f -m^ , U = {pi-pAf 
and has a perturbative expansion similar to Eq. ([3]), 



m 



In terms of the amplitudes the expansion coefficients in Eq. ([8]) may be expressed as 



(0)|fAf(0))=2^^^^(l-£) 



ut 



2— (r + u^- £r ) + (m) 



(6) 
(7) 

(8) 

(9) 
(10) 
(11) 
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where we have discarded powers in the heavy-quark mass mm Si . 

The expressions for Si^ have been presented e.g. in Refs. [8,9] and the loop-by-loop contri- 
bution (fW (i^lfW (^)) at NNLO in has been published in Ref. [31]. Both resuhs for and 
(|<7Vf (1) I (1)^ have been obtained in dimensional regularization and with the complete dependence 
on the heavy-quark mass. Here we provide for the first time the real part of (fAf |5Vf (2)) up to 
powers O (m) in the heavy-quark mass m. 

3 Method 

3.1 The massive amplitude from QCD factorization 

Let us briefly recall the key findings of Ref. [25] on how to calculate loop amplitudes with massive 
partons from purely massless amplitudes. Heuristically, the QCD factorization approach rests on 
the fact that a massive amplitude M for any given physical process shares essential properties in 
the small-mass limit with the corresponding massless amplitude M ("'=0). xhe latter one, M ('"=0)^ 
generally displays two types of singularities, soft and coUinear, related to the emission of gluons 
with vanishing energy and to coUinear parton radiation off massless hard partons, respectively. 
These appear explicitly as poles in £ in dimensional regularization after the usual ultraviolet renor- 
malization is performed. In the former case, the soft singularities remain in 5Vf ('"^ as single poles 
in £ while some of the coUinear singularities are now screened by the mass m of the heavy fields, 
which gives rise to a logarithmic dependence on m, see e.g. Ref. [32]. 

This structure of singularities for massless amplitudes has been clarified to all orders in pertur- 
bation theory [24] as all l/£ terms can be exponentiated, see also Ref. [33]. Similarly, all poles 
in £ and log(m) terms for amplitudes with massive partons also obey an all-order exponentiation 
with mostly the same anomalous dimensions as in the massless case [25]. Thus, in the small-mass 
limit the differences between a massless and a massive amplitude can be thought of as due to the 
difference in the infrared regularization schemes. QCD factorization provides a remarkably simple 
direct relation between Uv[ and M ('"=o) 

5ifH = Y[ X 5if(-=o). (12) 

(£ {all legs} 

The function Z^'"!*^) is process independent and depends only on the type of external parton, i.e. 
quarks and gluons in the case at hand. For external massive quarks Q it is defined as the ratio of 
the on-shell heavy-quark form factor and the massless on-shell one, both being known [34-36] to 
sufficient orders in tts and powers of £. An explicit expression for 

up to two loops is given in Ref. [25] (note the different normalization tts/ (47l) used there). The 
leading rir terms ~ {nM^Y for our process Eq. ([T]), gg QQ, can also be predicted based on the 
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Z,r = l + S4>fS)'4>«(a^). (14) 



above arguments. Keeping only terms quadratic in and/or = nf^ + rij one has up to two loops: 
where 

4? = ^'^/7>4? + ^K' (15) 

with z\^J^ ~ given in Ref. [25] (again note the different normalization tts/ (4;r) used there). z\^J^ 
is also equal to the O (tts) term in the gluon wave function renormalization constant Z3 in Eq. (|2l)) . 
The relation of Zj^^"'*^"* to Z3 is discussed after Eq. (|24l) below. To derive Eq. ([T?] ) we apply the 

definition for Z^"'!^) given in Ref. [25], i.e. evaluate the ratio of the gluon form factor with heavy- 
loop insertions and the pure massless gluon form factor [36,37]. The additional renormalization 
constant that enters the effective Hgg vertex (see e.g. Ref. [37] for details) cancels in the ratio and 
does not contribute to Exploiting the predictive power of the relation Eq. (fT2)) and applying 

it to the process Eq. ([I]) we get 

2Re(fAf(0)|5if(2))W = 2Re(5Vf(«)|;^(2)^(™=0) ^ (^2f^)^zf^y)^6,(,.=o) 

+ 2 (zg +Zf;/ +Z[^}z[;)) ^4,('n=o) + o{n,' X n,«) + 0{m) , (16) 

which assumes the hierarchy of scales <^ s,t,u , i.e. we neglect terms 0(m). Eq. ([T6l) predicts 
the complete real part of the squared amplitude {9V[ ^^^\!M ^'^^)^'"^ except (as indicated) for those 
terms, which are linear in rij^ (the number of heavy quarks) and, at the same time not proportional 
to ni (the number of light quarks). These two-loop contributions have been excluded explicitly 
also from the definition [25] of Z^™!^), as one needs additional process dependent terms for their 
description. 

The two-loop massless amplitudes Re (fW liW (2)^(ni=0) computed in Ref. [26]. We have 
checked that the finite remainders of the squared two-loop amplitudes obtained after the infrared 
subtraction procedure discussed in that reference agree with the corresponding terms constructed 
from the two-loop helicity amplitudes calculated in Ref. [27]. We have also found similar agree- 
ment between the finite remainders of the qq q'q' amplitudes we used in Ref. [22] that we 
extracted from Refs. [38, 39]. 



3.2 Direct computation of the massive amplitude 



The direct computation of the massive amplitude proceeds according to the same scheme as in 
our previous publication [22], which itself was an evolution of the methodology developed in 
Refs. [40-42]. In short, the complete amplitude is reduced to an expression containing only a 
small number of integrals with the help of the Laporta algorithm [43] . In a next step, Mellin-Bames 
(MB) representations [44, 45] of all these integrals are constructed, and analytically continued in 
the dimension of space-time with the help of the MB package [46] revealing the full singularity 
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structure. A subsequent asymptotic expansion in the mass parameter is done by closing contours 
and resumming the integrals, either with the help of XSummer [47], or the PSLQ algorithm [48]. 

We shall now concentrate on the differences with respect to our previous calculation. First 
of all, the number of master integrals is substantially larger, reaching 422, which adds a lot to 
the complexity of the computation. This is partly due to the fact, that the symmetry with respect 
to the exchange of the gluons generates the same topologies in the t- and w-channels, but more 
importantly because of completely new topologies, which come together with the more extended 
set of gluon interactions. Whereas in Ref. [22], it was possible to avoid the computation of the 
high-energy asymptotics of non-planar graphs, and still have a test of the factorization approach, 
we were not able to avoid them here. In fact, this additional complication is due to the single heavy- 
quark loop diagrams of Fig. [U which are explicitely removed from the factorization approach. The 
complete set of non-planar master integrals belonging to this class is depicted in Fig. |3l 




Figure 1 : Most complicated diagrams of the 
single heavy quark loop contribution. The 
thick lines are massive. 
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Figure 2: Definition of the momentum in the 
numerator of the non-planar integral consid- 
ered in the text, together with the labeling of 
the denominators. 




Figure 3: Non-planar master integrals cor- 
responding to the second diagram of Fig. [T] 
The numbers denote the multiplicities of the 
integrals of the given topology. All of the 
required MB representations can be derived 
from the MB representation of the integral 
Fig. [21 



The first problem that one has to face when dealing with non-planar integrals is the construction 
of MB representations. In the planar case, the iterative loop-by-loop integration has proved to be 
the most fruitful. On the other hand, the first non-planar double-box diagram ever computed [45], 
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with massless propagators and on- shell external legs, had its four-dimensional MB representation, 
l^P^' derived directly from the two-loop Feynman parameter representation. It seems, how- 
ever, that when masses are involved, the loop-by-loop representations are more compact, as seen 
for example in Ref. [49, 50]. However, any asymptotic expansion contains a so-called hard part, 
which is obtained by setting all the small parameters to zero, and which in our case would corre- 
spond to the massless graph. Following this line of thought, one can derive a representation for 
the massless on-shell graph by taking suitable residues in the result presented in Ref. [49]. Un- 
fortunately, one arrives at a six-fold representation, I'^p^'^^""^^ , in clear disadvantage with respect 
to I'^p'^' '^""^'^ . Interestingly, there is an even more severe problem inherent in the loop-by-loop 
approach. The leading pole derived in Ref. [45] reads (up to normalization factors irrelevant for 
this discussion) 



whereas the six-fold representation gives (with the same normalization) 



wr!=0, dim=6 
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This clear discrepancy is only explained when we look at the subleading pole from the six-fold 
representation, which contains a logarithmic singularity 

1.^-.-. ^ -4-^og{-s-t-u)+o (^) . (19) 

e^stu \^ J 

Obviously, the extension of the integral into the Euclidean domain, performed with the help of 
the u parameter, regularized part of the infrared singularity. The only way to obtain the correct 
result would be to first take the limit u —s — t, and only then £ — 0. This, however, is a highly 
non-trivial task, as is well known from studies aiming at the derivation of exact expressions of 
Feynman integrals in J-dimensions. 

In view of all the above arguments, we derived our MB representations directly from the two- 
loop Feynman parametric representation. In particular, for the scalar integral of Fig. |2l we have 
the following six-fold representation 



i=l 

xr(-z2)r(-z3)r(-z4)r(-z5)r(-z6)r(-2zi-z2-z3-2z4 + i)r(-£-z4)' 

xr(zi +Z2+Z3+Z4)T{-e + Zl +Z2+Z3+Z4-Z5 - l)r(-2£ + Zi - Zg - 2) 

X r( -2£ + zi+Z2-Z5-Z6-2)r{-2e + zi+Z3-Z5-Z6- 2) 

xr(z2 +Z6 + l)r(z3 +Z6 + l)r(z5 +Z6 + l)r(2£ -Zi +25 +Z6 + 3) 
x(^r(-2£-2z4)r(-3£-Z4-l)r(-2£+Zi+Z2+Z3-Z5-l)^) , (20) 

where the loop integration is done with the measure e^^^ / k / {in'' per loop. We defer the 
presentation of the full result of the expansion of this integral to the Appendix |Al Here we only 
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note, that the leading term of the expansion has a square root singularity, which is a feature of 
non-planar graphs, that does not occur in any of the planar integrals considered in this calculation. 
Clearly, the disappearance of this square root singularity is a simple test of the correctness of the 
calculation. 

The second problem that requires care is connected with the choice of the master integrals. In 
Fig.|3l we have not only shown the topologies, but also the multiplicities of the masters. The basic 
seven-liner needs as much as five different integrals. It is clear that we want to avoid coefficients 
containing poles in £ or m^, since the leading behavior of the integrals is difficult enough to deter- 
mine, and such poles would be synonymous of higher orders in the respective expansions. After 
inspection it turns out that one can take two integrals with second powers of the denominators into 
the set, but tensors rank one and two are unavoidable. Although it is possible to generate represen- 
tations for arbitrary tensors, one ends up with a huge number of four-fold integrals after expansion. 
Instead, one can introduce a new ficticious propagator that will have a negative power. For this we 
choose the square of the momentum that runs through the crossed box subloop, as shown in Fig.|2l 




Figure 4: A suitable choice of integrals that forms together with the underlying scalar integral the 
set of five master integrals of the seven line non-planar topology. The momentum W has been 
defined in Fig. [21 whereas a dot on a line denotes a squared propagator. 



The final set of seven-line non-planar master integrals is shown in Fig. HI Note that the dotted 
masters, i.e. those having higher powers of chosen propagators are particularly easy to calculate, 
because one only needs the singularities in l/m^. All of the representations can be derived from 
the following one 

xr(-Z3)r(-Z4)r(-Z5)r(-Z6)r(-Z7)r(-Z8)r(-Z2 + Z4+Z5)r(zi+Z2+Z3+Z6) 

X T{a -I- 2e - Zl -I- Z7 Z8 - 4)r(z2 + Z2,-ZA-Z5+Z2,+a \)T{z5 + zs + a^Wizi + Zi+a^) 

X r(-2zi - 2Z2 - 2Z3 + Z4 + Z5 - 2Z6 + <37)r(z2 - Z5 + «8)r(-£ - Z6 " + 2) 
X r(-£ - Z6 - aiA + 2)r(-2zi - Z2 - Z3 + 24 - 2z6 + ai%) 

X r(-£ + Zl + Z2 + Z3 - 24 + 26 - Z7 - «5678 + 2)r( -2£ + Zl - Z8 - «1234678 + 4) 
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xr(-2£ + Zl +Z5 -^7 -Z8 -'31245678 +4) 

X r(-2£ + Zl + Z2 + Z3 - 24 - 25 - 27 - 28 - « 1345678 + 4) 

X (r{-a -3£-Z6 + 6)r(-2zi -zi- 2z3 +24 - 2z6 + a78)r(-2£ - Ize - a 1234 + 4) 

X r( - 2£ + Z 1 + Z2 + 23 - 24 - 27 - ^2 1 35678 + 4) 

xr(-2£ + Zi+Z2+ 23 -24 -27-^245678 + 4)) , (21) 



where a = Y^f^i ai, as = L/es with S a subset of 1 , . . . , 8, and a,- with i= 1 , . . . , 7 are the powers of 
the denominators according to the labeling given in Fig.O whereas ag is the power of the additional 
denominator l/k^. Even though Eq. ((2T)) has two more integrations than Eq. (|20|) . the presence of 
the factor 1 /T{a'g) makes it necessary to perform first an analytic continuation in ag to a negative 
integer value, which effectively reduces the number of integration variables back to six. In fact, 
the above representation can be used to compute any of the master integrals from the full set of 
non-planars of Fig. [31 

To complete our exposition of the direct calculation of the amplitude, we have to note that the 
renormalization of the bare amplitude requires the on-shell wave function renormalization constant 
of the gluon, which is non- vanishing due to the presence of heavy-quark loops. We give it here as 
an expression exact in J-dimensions, 



Z3 



1 +a^nhTf 
+Cf 



2 

3£ 

4£3-7£- 1 



+ ( a^] n^Tf (nhTf ( 



-4£5 + 15£3+£2-ll£-3 



^£(4£3-8£2-£ + 2) J V2e2(4£4-4£3-13£2 + 7£ + 6) 
where is defined by the bare coupling constant and with 5e = (4;r)^exp(— £Ye), and 



(22) 



271 V rn^ 



/Y£r(l+£)5e. 



(23) 



In terms of the renormalized coupling and expanding in powers of £ through sufficient terms, 
the result for Z3 reads 

-3 - 



18 

+niTp 

+Ca 



4 4 



9£ 



2 
3 



A* 



2-, 



27 



9£2 9£^°^Vm2 



A* 



A* 



I 21 - 



35 
L36£2 



18£ 



5 

8£ 



1 

9 



13 
48 





15- 






T. 




137r2i 




(24) 


^ 216 \ 


}■ 



The corresponding expression for the wave function renormalization constant Z2 of a light quark 
state has been given in Ref. [22]. 
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Our result for Z3 in Eqs. (I22l) and (l24l) coincides to first order in ttg with Z^^j from Eq. (fT4l) . 
At second order in tts all known terms (i.e. those quadratic in the number of flavors) in the two 
constants are also a complete match (also in J-dimensions). In our direct evaluation of Z3 in 
Eq. (I22I) we observe gauge independence through two loops (within the class of covariant gauges 
employed in the calculation), which is consistent with the arguments given in Ref. [25] in favor of 

(mlO) 

the identification of the constant Zr , with Z3 evaluated in a physical gauge. 



4 Results 



We are now ready to present our result for gg QQ scattering for the interference of the two-loop 
and Born amplitude, 

2Re (fAf 1 5Vf (2) ) = (A^2 _ ^ /^3^ ^i^B + l-C+ -^D + N^n,Ei + N^n^Eh (25) 

\ N N 

2 2 

+njFi + nj^Fh + ^Gj + ^G, +Nnj''Hj +Nnin^Hjh + Nn,^^Hh + + ^//, + ^4 



which we have ordered according to the power of the number of colors A'^ and the numbers of 
light and = 1 heavy quarks with nj: = n^ + rij^ total flavors. 

The coefficients A, Ei, Hi, Hih, Hh, Ii, Iih, and //, have been computed with both of our methods. 
We have found agreement between the direct computation of the relevant Feynman diagrams in the 
small-mass expansion and the QCD factorization approach as given by the universal multiplicative 
relation (fT2)) . All other terms linear in n^, that includes Eh, Fh and G/, have been obtained by means 
of a direct calculations of the massive loop integrals as detailed above. The remaining coefficients 
B, C, D, Fi and G/ have been derived by application of the factorization formula. 

We choose x= —t/s as the only dimensionless kinematic variable in the problem and we keep 
the dependence on the renormalization scale, /j, explicit. We also introduce the following compact 
notation 



L„, = log ( , = log (j^^ , U = log (x) , Ly = log ( 1 - x) . (26) 



The different components now read 



A = — <i -8a;^ + 8a; + 



+ 16- 



65 

6x 



l—x X 
I54x^ I54x 



-8;c'^ + 8;c + 



l-x 



+ L, 



2 I m 



2x —Ix- 



l-x 



—^ + 2-1 
2(1-^;) 2x 



+ 4 +Lv 8x-'-16;c+12-- + 



16a;2- 16;c 
65 



l-x 
142 



6(1 -x) 



+ LmL.v 



7 4 4' 
16x2 -16x-^ + 16-- 



1-x 



-16x^ + 16x + 



4 4 
16 + - 

l—x X 



+ L„ 



-22x2 + 22x + Lv 4x- 



l-x 



+ 2 



10 



+lJ4x^-?,x + 6--]+ ^ 



1 



+U -16x-^ + 32x-24 + ' 



2(1— x) ~^ 2x ~^ 
7 7 

-; +36-- 

\—X X 



Ux^ -AAx + U-i -16x2 + 



1-x 



+ ( -4x2 + 12x-9 + - ]l? 



+ ( 18x^-41x + 31-- )U- 



104x2 104x 
^^ + - 



+ L; -4x^-4x + 



W 8x2 _ _ 



1-x 



+ 8 



+- 



+16x + 



2x2 2x 



1 2 1 

+ ■ 



3 3 6(1 -x) 3 6x 



+ LJ 18x2 + 5x- 



+ L2l 



1-X 



+ 8 



1-x 
161 



161 100 

18(l-x) ~T8l'^~9~ 



-4x2 + 4x + 



1 



1 



1 4 + 
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Notice that apart from the classic polylogarithms up to weight four, our results are also ex- 
pressed in terms of Nielsen polylogarithms 



ri log"-i();)logP(l-x3;) 



(43) 



5 Conclusions 



In the present article we have computed the two-loop virtual QCD corrections to the production of 
heavy-quarks for the gluon fusion process in the ultra-relativistic limit. Previously we had already 
obtained the corresponding results for quark-quark scattering [22]. Taken together, these results 
complete the two-loop radiative QCD corrections to heavy-quark production in hadron-hadron 
collisions in the limit when all kinematical invariants are large compared to the heavy-quark mass. 
The remaining channel with a gq (or gq) initial state involves at most one-loop corrections in a 



27 



consistent NNLO treatment. At one-loop the complete structure of the singularities as well as all 
large logarithms in the heavy-quark mass m can therefore be entirely treated with the methods of 
Ref. [32]. 

Our derivation relies on the combination of two completely different methods and we have 
ensured substantial overlap between them. In this way we have had mutual and highly non-trivial 
cross checks on our direct calculation of massive Feynman diagrams and on the QCD factorization 
approach [25]. In particular through the latter method we have had the possibility to relate to 
various massless [26-28, 35, 36] and massive results [34] available in the literature and we have 
found consistency. 

The results for gluon fusion in the present article (and Ref. [22] for quark-quark scattering) are 
of direct relevance whenever power corrections in the heavy-quark mass are negligible. This is 
certainly the case for hadro-production of bottom pairs over a large kinematical range at colliders 
and, to a lesser extent perhaps, for ff-production at LHC energies. Here, possible improvements 
would come from the systematic computation of power corrections in the heavy-quark mass, i.e. 
terms proportional to (m^/s)^ with ^ > 1 to improve the convergence of the small-mass expansion. 
This could be achieved, for instance, by extending the methods of Sec. |3]for the direct calculation 
of massive Feynman diagrams to higher powers in m. 

Finally, it is clear, that our result for ^^^^M ^^^) still has to be combined with the tree-level 
2^4, the one-loop 2 — > 3 as well as the square of the one-loop 2^2 processes (fW ^^^^M ^^^) in 
order to yield physical cross sections. Some of the matrix elements including the full mass depen- 
dence can be easily generated, others have become available in the literature only rather recently, 
see e.g. Refs. [31,51]. The combination of all these contributions enables the analytic cancella- 
tion of the remaining infrared divergences as well as the isolation of the initial state singularities. 
The latter will have to be absorbed into parton distribution functions of, say, the proton in order 
to match with a precise parton evolution at NNLO [52,53]. All these remaining steps are neces- 
sary prerequisites e.g. to the construction of numerical programs which then provide NNLO QCD 
estimates of observable scattering cross sections for heavy-quark hadro-production. 

A Mathematica file with our results can be obtained by downloading the source from the 
preprint server [hFtp: //ar Xi v. org [ The results are also available from the authors upon request. 
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A Asymptotics of the massive non-planar scalar integral 



Here, we give the leading high energy behavior of the non-planar scalar integral with a massive 
loop corresponding to the MB representation Eq. (I20l) 
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